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I. INTRODUCTION 



Gyrofluid equations in toroidal geometry were originally derived for small-amplitude dis- 
turbances to treat linear instabilities 111. Nonlinear terms were taken from an earlier version 

n 

derived in slab geometry [2] . Landau damping dissipation, resulting from phase mixing into 



Tom a previous closure dissipation 
% |5| . The original idea behind gy- 



arbitrarily small velocity space structures, was taken 
model [sj. Electromagnetic versions were also derived 
rofluid equations, derived in two-dimensional slab geometry, was for the density of a fluid of 
gyrocenters rather than particles, with a polarisation equation tying the gyrocenter densities 
together in a statement of overall charge neutrality ^\. Polarisation results from the part 
of the space density arising from the gyrophase angle dependent part of the distribution 
function, while the gyrocenter density reflects the part independent of the gyrophase angle. 
This was already a result of gyrokinetic theory [7|, which was originally formulated as a 
finite Larmor radius (FLR) correction to the drift kinetic equation |8| . The gyrofluid model 
rests on the gyrokinetic model, much in the same fashion that the more familiar fluid model 
rests on the Vlasov (or Boltzmann) equation for particles 

The gyrokinetic theory itself, however, underwent a serious advance resulting from the 



use of Lie transformation techniques applied to the drift kinetic Lagrangian |llHl3|. The 
Lagrangian was to be transformed to gyrocenter coordinates and only then the gyrokinetic 
equation was to be found in terms of the resulting Euler-Lagrange equations describing the 
motion of individual gyrocenters. The gyrofluid model, on the other hand, continued to be 
based upon the older gyrokinetic formulation, in which it was not always transparent what 
one had to keep in order to maintain consistency at a given ordering. Indeed, the original 
toroidal model had to be modified before it had a clearly consistent energy conservation 



theorem jl^ . Il5| . The gyrokinetic theory was re-cast as a Lagrangian field theory 
and nonlinear nonlocal gyrofluid equations have also been formed in a similar manner [l8 



The delta-f form of the gyrokinetic theory also has its energy theorem 20|, whose con 



m 

19|. 



served quantity as a functional quadratic in all the dependent variables is better thought 



of as a free energy, similar to that of fluid equations |21|-|24|. This has been shown to be 
related to en trop y [25| . Free energy for fluid equations has also been shown to be related to 
this entropy [26]. 

The original gyrofluid derivation for toroidal geometry led to some inconsistencies which 
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caused fluctuation free energy not to be conserved. These were repaired by a construction 
method which can be thought of as somewhat arbitrary and/or artiflcial is!]. The free 
energy was determined by analogy to the fluid equations, leaving it open to the argument 
that it might not be fundamental. The key insight in the meantime was the connection 
made between the models, showing the gyrofluid moment variables can be cast in terms of an 
Hermite polynomial representation of the gyrokinetic distribution function. The polynomials 
are functions of the velocity space coordinates, and the coefficients of these are the gyrofluid 
moment variables. Given the gyrokinetic delta-f free energy functional, the gyrofluid one 



can be derived using this Hermite expansion. This was shown in Ref. 
simplifled gyrofluid model without FLR corrections. 



26( 1 in terms of a 



Herein, we systematise this procedure for the six-moment gyrofluid model of Refs. |2|, |4| 



which was corrected for energy conservation in Ref. 



15| . But before that we need to show 



the relation to the only really flrst-principle model in the hierarchy: the total-f gyrokinetic 



equation obtained by Lie transforms 12 



13| . We start with that one, simplifying it mildly 



to obtain a computationally tractable form. Then we make the delta-f approximations in 
such a way as to kee p energ etic consistency intact. The total-f energy obtainable from 



the Noether theorem 



16|, 



17| is replaced by the delta-f free energy functional referred to 



above. This delta-f equation and its free energy theorem then launch the derivation of 
the gyrofluid equations in local form (fully nonlinear equations, but constant background 
parameters, dependent on drift ordering). At the kinetic level, the gyrocenter part of the 
charge density in the polarisation equation is the same in the total-f and delta-f versions, 
and arises from the interaction term in the total-f Lagrangian. The FLR closure for the 
gyrofluid model is applied only once, and it is here. Ultimately, the role of polarisation is 
the same in all three models, and taking that into account allows energetic consistency to 
be maintained. Speciflcally, the same rules must be applied both to polarisation and to the 
moment variable equations 15|. The polarisation part of the charge density is the same 
in all models which use linearised polarisation, including the delta-f and gyrofluid levels. 
This obtains the polarisation equation, which at all levels can be used to re-cast the ExB 
kinetic energy in terms of gyrocenter charge potential energy. The part of the gyrofluid free 
energy due strictly to the moment variables follows from insertion of the Hermite polynomial 
representation into the delta-f form. The rest of the derivation is a mere matter of consistent 
application of the moments to the delta-f gyrokinetic equations. This part of the method is 
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similar to that of Refs 



2|, but the results of Ref. jl5| concerning the conservative part 
of the equations and the free energy are recovered without further pitfall. Dissipation is 
inserted manually as before. 

Previous versions of gyrofluid equations have applied the FLR closure only to the electro- 
static potential, leaving the electromagnetic response along. This is usually justified since 
the turbulence spectrum only extends down in scale as far as the ion gyroradius. Hence, 
FLR electron effects are small by the mass ratio. Moreover, ion dynamical contributions to 
the Ohm's Law are also finite mass ratio corrections. However, fluid-type equations such 
as these are also useful in the study of coUisionless reconnection, which includes sufficiently 
small-scale phenomena that electron FLR effects may enter 27l430l|. Herein, the FLR treat- 
ment is extended to the magnetic potential at the same level of sophistication as to the 
electrostatic potential. The parallel velocity and perp/parallel heat flux moments are mixed 
in the magnetic flutter disturbance and induction physics in the same way as the density and 
perpendicular temperature moments in the ExB advection and polarisation physics. This 
extension of the gyrofluid model represents the new result of this work. 

II. THE TOTAL-F AND DELTA-F FORMS OF THE GYROKINETIC MODEL 

The total-f Vlasov equation is a Hamiltonian bracket equation for the total distribution 
function /, a dependent variable over phase space coordinates {x, z, w} with the space part 
x given by field aligned coordinates {x,y,s}. Field aligning refers to a single nonvanishing 
contravariant component of the background magnetic field, in this case B'^, which identifies 
s as the parallel coordinate. The velocity space part is given by the parallel velocity z 
and magnetic moment w, with / independent of the gyrophase angle 6. As a starting 
point we use the result of the Lie transform theory version of the gyrokinetic expansion. 



specifically a version of the forms given in Eqs. (16,17,22-24) of Ref. 13[, mildly simplified 
for computational use. The particle Lagrangian is given by 

Lp = -A* ■ — + w r-H (1) 

^ c dt e dt ^ ^ 



where 



A* = A + mz^h (2) 
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where A and b and B are the potential, unit vector, and field strength of the background 
magnetic field. The Hamiltonian H is given by 



H = m^ + wB + eJo(j) - ^g^iM^I)'') " 



(3) 



where is the electrostatic potential, Jq is the gyroaveraging operator, and the contribution 
which is quadratic in is referred to as the gyroscreening potential. Formally, Jq has the form 
in wavenumber space of multiplication of Fourier coefficients by the zeroth Bessel function 
Joik^PLj, where pl is the particle gyroradius given by pi = v_}^/VL, with gyrofrequency 
Vt = \eB/mc\, or in terms of the coordinates by p\ = 2wB/{mQ'^). The particle equations 
of motion are found from the Euler-Lagrange equations resulting from Lp, 



Bt 



zW 



-F-VH 



dt eB 
Drift tensor notation is used, with definitions 



mB\ 



. dz 
'dt 



-W ■ VH 



(4) 



B* = VxA* 



Bl 



b B* 



F = VA - (VA)^ 



It follows that 



F = e B 



Vxb 



c F 

B* = B - mzV 

e B 



(5) 
(6) 



(7) 



where e is the rank-three Levi-Civita pseudotensor. 

The polarisation equation is found by extending the Lagrangian to obtain one for the 
entire particles/field system using the methods of Refs. [l6|, Il7j . 



L = Y.fdAf 
sp 



e . ^ dx mcdO 
-A* ■ — +w - H 

c dt e dt 



(8) 



where the sum is over particle species, / is the distribution function for particles of each 
species, Lp multiplies / expressed in phase space coordinates, H is given in Eq. ([3]), and the 
phase space, velocity space, and configuration space integration domains are given respec- 
tively by 



dK 



dV / dW 



dV = d-^x 



dW = 2nm'^dzdwB, 



(9) 



hence identifying Bf! as the variable part of the velocity space volume element. 
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The Euler-Lagrange equations for are found by varying the Lagrangian with respect 
to 0, yielding an integral over dV of 50(x) times a coefficient, which is required to vanish. 
This produces 

Y^JdW [eJof + Jo [MM) - (JoM) <P] = (10) 



sp 

where Ai is the polarisability given by 



B ow 



(11) 



It is customary in most treatments to linearise the polarisation term involving Ai. This 
corresponds to replacing / by a Maxwellian background F'^^ in the term in the system 
Lagrangian which is due to the gyroscreening potential. Since the dependence of upon 
w is proportional to exp{—wB/T), the polarisability is replaced by 



M 



e'F^'/T 



The Hamiltonian is replaced by 

H 

and the polarisation equation is replaced by 



m— + wB + e Jo<; 



E 

sp 



e Jo/ + Jo 



T 



-Jn 



Jn 



T 







(12) 



(13) 



(14) 



where now only the source term depends on the variable /. The facts that H is linear in (j) and 
polarisation is linear in / lead to the fact that this particular set of terms remains intact when 
the delta-f approximations are taken. Ultimately, mutual field/particle energy conservation 
depends on the properties of these terms, and the gyrofluid energy considerations descend 
directly from the total-f gyrokinetic ones. 

Most treatments neglect the action of Jo upon F^ , leaving the familiar form with 
(e^/T)F^^(jQ — 1)0, whose velocity space integration produces (ne^/T)(ro — 1)0, where 
Tq is an operator whose form in wavenumber space is multiplication of Fourier coefficients 
by Io{h)e~^ with argument h = k'^pf evaluated with the thermal gyroradius pi, in turn given 
by pj = MTc^ /e^B"^ for species i [7|. In recent total-f gyrokinetic computations, this op- 
erator is further simplified by using V ■ {nMc^ / B'^)V(f), its low-A;_|_ limit 3l|-|33|. In these 
forms n and T are background constants or, at most, radial (across magnetic flux surfaces) 
profiles for the species density and temperature. 



When polarisation is linearised, the corresponding term in the Lagrangian becomes a 
field energy term, with / replaced by F^^ . Hence the second order field piece (gyroscreening 
potential) is moved from H into the field energy term in the Lagrangian. Consistency 
is maintained only when these steps are taken together: if polarisation is linear then the 
gyroscreening potential must be absent in H, or vice versa 16 1. 

The equation for / is now determined by the particle equations of motion with the 
linearised if, with the actual form remaining intact. It is given by 

Bf-l + Vii ■ -| ■ V/ + -B* ■ f ^V/ - ^^Vh] = (15) 

" ot e B m \oz az ) 

Eqs. ( !T3|14fl5p form the electrostatic version of the total-f gyrokinetic model which is our 
starting point for the rest of the derivation. First, the simplest extension to include shear 
Alfven dynamics is given, then the delta-f model and its energy theorem are presented, and 
then the gyrofiuid model descends from that. 



A. Electromagnetic extension 



The discussion of the electromagnetic version of the model is kept brief because the 
point is not to discuss various representations but to establish what is needed to proceed to 
the delta-f formulation. This version retains z as the parallel velocity coordinate, not the 
canonical parallel momentum (cf. Ref. 134]). The particle Lagrangian Lp retains the form in 
Eq. ([1]) with A* now expanded to include the parallel magnetic potential, Ay, only, so that 

(16) 



A* = A + ( JoAii + mz- 1 b 



The main complications are that the part multiplying b is not only spatially variable but 
contains a time dependent field variable. The particle equations of motion become 



Bl— = zW-—Y-^E 

eB 



mB;^ = -B*-VH--Br,Jo^ (17) 



II dt eB ' " " dt 

with B* now containing magnetic nonlinearities involving V(Jo^||)- The overall Lagrangian 
including the fields becomes 



L = J2fdAf 

sp 



(ix mc dO 
— + w — 

dt e dt 









H 


-fdV^ 


V±A|| 




J Stt 





(18) 



with the last term representing the magnetic energy, and with A\\ also occurring in A*. 
Variation of this with respect to is unchanged, and produces Eq. (|T0|) or its simplifications 



such as Eq. (fT^ . Variation with respect to produces the induction equation (gyrokinetic 
version of Ampere's law), 

Ajr f 

VlA|, = dW[eJozf] (19) 

c gp J 

where the term on the right gives the current, and without using the canonical parallel 
momentum as a coordinate there is no complication with the skin depth. The latter does 
however remain implicitly present as the electromagnetic gyrokinetic Vlasov equation for 
this Lagrangian, 

- + Vff . H| . V/ + 1b- . f ^V/ - §^Vff) = (20) 

" at mc " oz ot en m \oz oz J 

now includes the time dependent induction effect through the explicit time derivative of A». 



B. The delta-f approximations 

The delta-f Vlasov equation results from the expansion / = F^^ + 5f , with 5f a correction 
to a Maxwellian background whose sole spatial dependence is through the magnetic field 
strength i?, and application of the delta-f ordering. Prior to the application of strict delta- 
f ordering, the parallel dynamics is linearised, so that the only nonlinear terms are those 
resulting from VH ■ F ■ V or B* ■ V, arising from the perturbed contributions of iJ, A\\ and 
/. With pijL ordered small (with L any background scale), B*^ is replaced by B. Further, 
the curvature is set such that b • Vb = log 5. This yields 

b| + V*. . v( V) + !!!£!±i^£ V log B 4f . v;. 

ot B e B 

+ l^B.V/.-l^B.V/7o = (21) 

m oz m oz 

where the auxiliary variables are 

g = 5f + —e-JoA\\ h = 5f + —eJo<j) = Jo - -^11 j (22) 

These are the inductive response, nonadiabatic response, and gyrokinetic potential, respec- 
tively. The zeroth order Hamiltonian and the Maxwellian are given by 

Ho = m— + wB F^' = n(27rT/m)-3/2 exp(-/Jo/T) (23) 

with m and e and n and T the species background parameters. 
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The terms in Eq. (l2Ti) are referred to as nonlinear advection, magnetic drift, and parallel 
trapping and streaming, respectively. Nonlinear advection is by the perturbed Hamiltonian, 
which includes ExB advection and the magnetic flutter nonlinearities, through and A\\, 
respectively. The grad-B and curvature drifts are combined, and departures from an inverse 
major radius dependence of B are neglected. The parallel dynamics includes both streaming 
{dHo/dz) and magnetic trapping (VHq = wBVlogB). The parallel dynamics is linearised, 
with both static and inductive pieces of E\\ appearing with {mz/T)F'^^ , leading to the 
expressions h and g, respectively. 

Now the strict form of delta-f ordering is applied. This refers to the split between the 
perpendicular and parallel coordinate directions following from the k\\ ^ kj_ ordering. A 
field aligned coordinate system is assumed, with x and y the perpendicular coordinates 
respectively following the radial and electron drift directions, and the coordinate s following 



the parallel direction. This is done for a general tokamak geometry following Refs. [35l-l37|. 
Axisymmetry renders d/dy = for the background magnetic field B. The unperturbed 
parallel gradient is given by 

d Odd 
B-V = B'— hence ~ TT > TT (24) 

OS ox Oy OS 

Under this ordering we have B = B{s), and hence B and F*^ are independent of x and 
y. However, log-B in the magnetic drift terms must be re-formed in order to keep both the 
interchange (dlogB/dx) and geodesic [dlog B /ds) curvature contributions. The derivatives 
of log B are evaluated, and then the resulting form is restricted to have dependence upon s 
only as with the rest of the geometry. This introduces the curvature operator /C, which is 
expressed in terms of a vector-gradient contraction, 

+ (25) 



where the components are given by 



rF 

)Ci^,y} = -^logB'- — -V{x,y} (26) 



It is necessary for free energy conservation that /C be a pure divergence. This is guaranteed 
by maintaining 

d d 

— /C" + —icy = (27) 
ox oy 
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which is trivially satisfied by having the /C" components depend on s only. Correspondence 
to the linear forms in Refs. is to identify /C with —2iujd, where Ud is the toroidal drift 
frequency. 

After this, the derivatives in the nonlinear advection term are restricted to {x, y} only, 
and the trapping/streaming terms involve derivatives in {s,^} only. The resulting delta-f 
ordered gyrokinetic equation is 







where the brackets are given by 



[f,9] 



ab 



df dg df dg 



(28) 



(29) 



dx°- dx^ dx'' dx°- 

with a and h denoting two of the phase space coordinates, /C is the magnetic curvature 
operator. The coefficients c, m, e, and Bq are constants; x, s, and w are coordinates; 
and -B, the metric coefficients contained in Jq, and the components of /C are functions of 
s only. The Maxwellian F*^ is a function of 5, z, and only. Field aligned Hamada 



coordinates constructed as in Refs. [361, 1371 are assumed. Hence, B^ and F^^ / B"^ are flux 
functions, which in the local limit reduce to constants. Following this is the use of the 
constant i?o, the average value of 5 on a flux surface. 

The fleld variables (f) and A\\ are functions of x and are determined in the total-f theory by 
the polarisation and induction equations, respectively. In the total-f theory these equations 
are found by variation of the total action with respect to the fleld variables. The delta-f 
forms are the corresponding linearised versions. The polarisation equation is 



E 

sp 



pM 

eJ^g + e^ — {Jl-l)(l) 







resulting from Eq. ( fT4l) . The induction equation is 

47r 



sp 



2 TPM 

ezJog - —z 

c T 







(30) 



(31) 



resulting from Eq. f|T9|) . In each case the delta-f form of the velocity space integral / 
has reduced to 

dyV = 2Txm~^ dz dw B (32) 

The source terms in these equations follow directly from the lowest order interaction La- 
grangian terms. 



dA {5f){-eM) 



I dK{5f)'-{eJoA\{) 



(33) 
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respectively, where we can replace / by Sf because these terms are linear in /. The Hermitian 
property of Jo and its commutation with z plays the central role. This is why the way Jq 
mixes moments through its w-dependence must be established the same way for the field 
variable derivatives in the Vlasov equation and the moment source terms in the polarisation 
and induction equations. The interchangeable use of g or 6f in Eq. fl30p follows from the 
antisymmetry of zJq in velocity space. 



C. Energy and free energy in the gyrokinetic models 



The conserved energy S for the total-f model is given by application of the Noether 
theorem to its Lagrangian 161,1121 • With H and L given in Eqs. (!3]|8|l . the energy is merely 
the phase space integral of H over /, 



S = J2jdAHf = Y.JdAlHo + eJo0 - ^^[^o(0') - W?]] f 

o [J o [J \ J 



(34) 



(35) 



sp " sp 

The polarisation equation (Eq. [10]) may be used to recast this in terms of 

M 

I dA ^ Hof + 

sp 

where the combination involving is the generalised ExB kinetic energy (in the low-k± limit 
it reduces to the familiar form in terms of v^, i.e., the usual ExB velocity squared). Note 
the use of the Hermitian property of Jq as well as integration by parts of d/dw. 

In the form with linearised polarisation, the energy is equivalently given by Eq. ( !34l) with 
/ — i- F*^ in the gyroscreening term. Again the polarisation equation (here, Eq. UM is used 
to eliminate the /eJo0 term, to find 

r r P*^ /A2" 

S = J: dA Hof + e^^{l-J^)^ 
sp L ^ . 

with Ai = {e^ /T)F^ inserted explicitly and the action of Jq upon /T neglected. The 
notation Jq0^ is shorthand for (Jo0)^- The velocity space integral over F^^ can be done 
explicitly to find 



(36) 



8 = Y.jdA 
sp 



/fo/ + n-(i-ro)|- 



(37) 



where again Fq^ is shorthand for (piVQcj)). 
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This field energy term is the same in the delta-f model. The interaction Lagrangian 
term (implicit in the delta-f model) is also the same, due to the linearity in both / and 
. This fact was used to motivate the discussion of free energy in delta-f models generally 



20, |25 



26| . Here, however, we have an even easier path, noticing that every term in Eq. 



is a bracket with h, i.e., a first order derivative on h and also a total divergence. It 
is trivial to multiply Eq. ( 128|) by h and any other quantity which commutes with all the 
brackets to obtain conserved quantities, called Casimirs of the system. The one identified 
with free energy is quadratic in h and must have dimensions of nT integrated over space. 
Hence the factor {T/F^'^) which works since depends on {szw} only and commutes with 
the 2;s-bracket, and T is constant. The contribution quadratic in 6f is much like a small- 
amplitude form of the thermodynamic entropy, while the contributions due to (p and A\\ 
become equivalent to the ExB and magnetic field energies present in the total-f model, once 
the polarisation and induction equations are used. Details and underlying considerations 
are in Ref. 251]. 

The energy theorem for these equations fl28ll30|31l) is given by 



'dt 



d_ 

dt 



j:[dA 

sp 



pM 2 



(38) 



or equivalently. 



dA 



?M 



pM 2 



T 



d 



dV — 

Stt 



Vi A, 



(39) 



reflecting mutual conservation of thermal and kinetic free energy, the E-cross-B energy, and 
the magnetic fluctuation energy, each given by 

T {5ff 



E 

sp 



dk 



pM 2 
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E 

sp 

M = 



r-M 



dAe^ 



dV 



T 



V, All 



(40) 

(41) 
(42) 



respectively. 

The need to use the field equations to bring the factor of h under the d/dt in Eq. ( I38l) is 
the reason the choice of second-order Casimir to identify as the total free energy is unique. 
Any integrable function of h multiplying Eq. (I28l) will produce a vanishing phase space 
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integral. However, the operation through Eqs. fl30|3ip will result in a single time derivative 
only if the multiplier is (T/F*^)/i. This yields the identification of {T/F^^)hg/2 as the free 
energy density and therefore S given by Eq. (15^ as the total free energy. 

III. THE GYROFLUID MOMENT SET AND FREE ENERGY 

The basic prescription of a gyrofluid model is the set of moments kept as dependent 
variables. The simplest version uses the density only [g^. One level up from that is a 
three dimensional version usin g de nsities and parallel velocities, equivalent in scope to the 
familiar four field fluid models jlj]. With the temperature gradient of either species setting 
the basic dynamical character, at least four moments, one each for the density n, parallel 
velocity u\\, and parallel and perpendicular temperature Ty and Tj_, are needed |2|. These 
temperatures are given by the parallel and perpendicular energy moments (over mz'^/2 and 
wB, respectively) divided by the density. The perp/parallel separation in the temperatures is 
made necessary by the underlying dynamics: only Ty enters parallel (hence Alfven) dynamics, 
and only T± is involved in polarisation (through the wB dependence of Jo). In the nonlinear 
dynamics there is moment mixing between n and T± but not with Ty , so even in the absence 
of curvature and grad-B drifts the responses of Tj_ and T\\ to the rest of the dynamics are 
different. Finally, magnetised plasma turbulence takes place at time scales for which the 
parallel sound wave transit frequency Cg/qR is very slow. Both the parallel viscosity and the 
perp/parallel components of the parallel heat flux have time-dependent responses to velocity 
and temperature gradients. In a gyrofluid model the parallel viscosity is proportional to the 
difference Ty —T±, which is already taken care of as Ty and T± have their own time dependent 
equations. This means that the perp/parallel components of the parallel heat flux (gyy and 
gj_y, moments over and 2:^5, respectively) also require their own time dependent 



equations. The six moment models 



15| are the result. 



Here and below, the tilde symbols are used to distinguish the dependent variables from the 
parameters. There is one set of moment variables per species, with all species contributing 
to polarisation and induction and the field energy pieces. 
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The moment variable list is given by 



n 



J dW (Sf) 
nu\\ = [ dWz{5f) 
nTii = I dW [mz^ - T) {6f) 



dW {wB - T) (Sf) 



gil= JdWimz'-3T)^{Sf) 
q^ = j dW {wB - T)z {6f) 



(43) 



where here and below q\\ and q± are used as shorthand for q\\\\ and q±\\, respectively. Heat 
fluxes perpendicular to B are given by the curvature and grad-B drifts, as combinations of 
n, T||, and T±, and therefore are not written explicitly. 

Use of a finite set of moment variables to represent 6f implies a representation of 6f 
in terms of a finite-degree polynomial in velocity space with those same moment variables 
as coefficients j26|. The Hermite polynomial decomposition of the distribution function in 
terms of the set of moment variables in Eqs. ( I43l) is given by 



Sf 



n u\\ z IT\\ \ T± 



wB 



1 



+ 



z 



-3 



z 



+ 



qj_ ( wB 



- 1 



z 
V 



(44) 



3nTV \ J V nTV \mV^ 

where V is the species thermal velocity given by V"^ = T/m. The coefficients are chosen 
for orthogonality and to recover the above definitions of the moment variables. Insertion 
of this form into the portion of the delta-f free energy in Eq. fl39l) dependent on 6f yields 
the gyrofluid free energy. The part dependent on the state variables, n or Ty or Tj_, is the 
thermal free energy, and the part dependent on the flux variables, u\\ or gy or q±, is the 
kinetic free energy. The spatial density of the thermal free energy is 

2 / ~ X 2" 

(45) 



1 /T|| 

2 \¥ 



T 



and the density of the kinetic free energy is 



U = — 
2 



3 \nTV 



+ 



nTV J 



(46) 



Together, Ut + represent the delta-f thermal free energy given by £f in Eq. (HOj) . and are 
derived directly from it. 
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The field energy pieces are the same as in the delta-f gyrokinetic version. Evaluating the 
integral over (1 - Jo^)F^ in Eq. (gl]), the gyrofiuid E -cross-B energy density is found. The 
magnetic energy density carries over directly from Eq. f l42p . These are given by 



E 

sp 



ne 



T 



ViAi 



(47) 



respectively. The operator Fq reflects gyroscreening. As noted above, its form in wavenumber 
space is the function 

ro(6) = h{h)e-^ with h = klp^ (48) 

where p = V/Q is the thermal gyroradius. Recall that in the local model Tq involves 
derivatives with respect to x and y only, while B = B{s) so that the Hermitian property is 
maintained. If wavenumber space is unavailable, the Fade approximant Fo(6) = (1 + b)~^ 
is used. Eqs. ( HSMTj) are the same as in Ref. [l4], but now they have a firm derivation in 
terms of the delta-f gyrokinetic version, applying the procedure of Ref. 261] to the result in 
Eq. ([39]). 



IV. CLOSURE RULES FOR GYROAVERAGING 

The evaluation of Jg in polarisation was trivial because it is integrated over F'^^ and the 
integral of Jo(6)^e~^ is well known. However, in the source term in polarisation (Eq. I30|) only 
one factor of Jq appears. A closure approximation for / dW JqF'^^ is needed. The one used 
previously was decided from the properties of resulting linear eigenfunctions . However 
it is possible to evaluate this directly. Taking the density moment over /i, we find 

r r F*^ 

J dWh = n + J dW —eM (49) 

In the second term, (f) is dependent on space only, and the fact that under the delta-f model 
F^^ commutes with perpendicular spatial derivatives may be used to do the velocity space 
integral separately. This moment defines the basic gyroaveraging operator, 

Ti = - [ dWF^^Jo (50) 
n J 

which acts the same way on any spatially dependent moment or field variable. It is merely 
a special function, as we find by inserting the form of F^^ , doing the integral over z, and 
defining x = wB/T, 

poo 

Fi(6) = / rfxe-Vo(v^) (51) 

JO 
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This form was used by Knorr ei a/ [g^. Conforming to the practice of Beer a/ however, 
we keep to the definition 

T,ib)^Tl/\b) (52) 
and hence the operation of Fi on a spatial variable is like that of Fq in wavenumber space, 

1/2 y 

with Fg (6) is used rather than Fo(6). If wavenumber space is unavailable, the Fade ap- 
proximant Fo(6) = (1 + h/2)~^ is used. The label Fi is used for generality: any Hermitian 
operator is admissible if consistency is the only requirement, so that if F*^ is non-Maxwellian 
then some form other than Fq''^ is chosen. These considerations are explained and justified 
in Ref. Q. 

In the polarisation equation the source term is / dW eJ^^df), which contains moments 
over both unity and wB times F*^. The corresponding closure approximation for this is 
found by inserting Eqs. (jHj) into the factor of eJo{6f) in Eq. (!30l) to obtain 



/ 



dWeJoiSf) = ne (v,^ + Fs^ j (53) 



where F2 is given by 



r2 = / dW^^^^^F'^Jo = / dWF^'Jo (54) 



following both appearances of T in the factors in which remain after the integral over z 
is done (cf. Eq. [23l) . It follows that whichever approximation is taken for Fi we always have 

where the latter form is the one to use in wavenumber space. 

We refer to Fi and F2 as the first and second gyroaveraging operators, or alternatively 
as the gyroaverging operator and its first FLR correction, respectively. For a six-moment 
model, this is as far as the moment hierarchy goes. If the moment over {wBY is a dependent 
variable in an extended version of the theory, then there is a further third gyroaveraging 
operator involving second derivatives of Fi with respect to b, and so forth. Strictly speaking, 
the {wBY moment will result in a form determined by d'^Ti/ d(\oghY . However, in this case 
we have to replace this 4th moment with a form mandated by energy conservation. 

Energy conservation in the delta-f Vlasov equation works with h appearing under all 
the derivatives and free energy evolution determined by multiplication of the equation by h 
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and integrating over phase space with sum over species (Eqs. [28]|38|39p . The corresponding 
form of this is that the moments of h appear under the derivatives in the gyrofluid moment 
equations. These moments are merely those over 5f together with contributions due to Jo0. 
Only two of these contributions are nonzero, as the others vanish due to odd symmetry or 
orthogonality in the integrals over z. The surviving ones are 



n + n-Tic 



wB - T , 

dW h 

T 



n— + n^V2^ 



(56) 



Hence in the linear terms in the moment equations, derivatives of n + n(e/T)ri0 must 
appear in that combination, and the same holds for Tj_ + er20. This was the procedure by 
which Ref. 15| repaired gyrofluid energy conservation. The higher moments (4th and 5th) 



over h which occur in the equations are then determined by the same ones over 5/, which 
are straightforward. Hence for the (wB)"^ moment we have 

(57) 



dW{wByh = J dW{wBy{6f) + ^ (^j dWwBwBF^^Jf, 
The first piece is found by straightforward evaluation 

J dW {wBf {6f) = 2pT P + 2 



T 



(58) 



Hence this moment over h must have the form found by combining Fi^ with n and with 



- + ;^ri^ 

n I 



J dW {wBf h = 2pT 
Subtracting these two we find the required form, 

Y dW wB wB F^^^Jo^ (j) = 2pT^ {Ti + 2T2 



e 
T 



(59) 



(60) 



and hence the requirement on the operators 



- dW 

n J \ I / 



wBy 



) F^Vo = 2 (Fi + 2r2) 



and also 



dW 



wB-T 
T 



F^-'Jo = (ri + 2r2) 



(61) 



(62) 



This replaces the form which would be found from two derivatives applied to Fi. But that 
third form is only viable in the case the 4th moments of 5f are retained as dynamical 
variables. In the six-moment model however the last equation above replaces that, and the 
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moment hierarchy is closed. The closure of higher moments over Sf is given by its form 
in Eq. (jS]), and then the higher moments over Jq are given by the requirement that the 
combinations in Eq. fl56p always appear intact in the moments over h. 

Equivalently, we list the two approximations involved in the gyrofiuid model: (1) the 
choice of moments to keep as dynamical variables, and (2) the form taken for Fi. The rest 
of the model then follows by simple evaluation constrained by energy conservation (keeping 
the pieces of moments over h together). 

V. GYROFLUID POLARISATION AND RELATION TO ENERGY CONSERVA- 
TION 

We now return to the gyrofiuid polarisation equation. Insertion of 6f from Eq. (jH]) into 
Eq. (!30|) and using the definitions of Fi and F2 in Eqs. fl50f54p . we have 



As in the delta-f gyrokinetic polarisation equation (Eq. I30l) . this is a statement of strict 
quasineutrality, with each species charge density given by the gyrocenter part (moment 
variables) and the polarisation part (due to the electrostatic potential). 

In the total-f or delta-f gyrokinetic models, we can recover the interaction Lagrangian in 
Eq. (!33|) by multiplying the polarisation equation by and integrating over space (essentially 
un-doing the steps by which the polarisation equation is derived in the first place). Under 
linearised polarisation the interaction Lagrangian is the same in both the total- and delta-f 
models, since the term is linear in both / and 0. We can do the same in the gyrofiuid 
version, either operation by / dV (px or simply by inserting 6f from Eq. (1441) into Eq. fl33|) . 
In all these forms, the Hermicity of Jq and hence Fi and F2 has a central role. 

For the result is 



sp L V JL /J 

whose variation with respect to recovers the gyrocenter source terms in Eq. fl63|) . Applying 
the Hermitian property of the F's, this is equivalent to 




(63) 




(64) 




sp 



(65) 
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where 



<Pg = t,{4>) 



n 



dW^^^^^F^'iM) (66) 



are defined as the first and second gyroaveraged potentials, respectively. This is the under- 
lying reason (pc is associated with n and with Tj_ and, ultimately, why they must appear 
together in derivatives representing energy transfer processes between the various equations. 

The corresponding forms for Ay use the delta-f induction equation (Eq. EI]) and the part 
of Eq. ([33D due to Au, with the result 



^int,j4|i 



sp •' 



n-A» riui, + r2 



T 



(67) 



Variation of this with respect to A\\ recovers the source term in the gyrofiuid Ampere's law, 

= (68) 



^ sp 



ne [Tiuii +T2Y 



We can also find this by inserting the form for Sf in Eq. ( 1^1) into the gyrokinetic Ampere's 
law (Eq. I3T1) and evaluating the velocity space integrals. Using the Hermitian property of 
the r's, the electromagnetic interaction Lagrangian can also be re-cast as 



sp 



n- [u\\Ag + — Xg 



(69) 



where 



Ag = Ti{A\\) = - / dWF^{J^A\\ 



xg = r2(^ii 



dW^^^^-^F^iJoA 
nl 



(70) 



are defined as the first and second gyroaveraged magnetic potentials, respectively. This 
leads to the appearance of u\\ and Ac together, and q± and xg, under d/dt in the gyrofiuid 
moment equations. 

All of this is closely related to conservation of free energy in terms of the functionals in 
Eqs. fl38|39l) for the delta-f gyrokinetic model and Eqs. (l45] - H7|) for the gyrofiuid case. The 
field energy components in Eqs. fH7|) can be re-cast using the polarisation and induction 
equations in Eqs. fl63f68p as 



e {(t)Gn + nnc'^ 



(71) 
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and 



M 



E 

sp 



dV 



^ {nAcUw +XGY 



(72) 



respectively. The time derivatives follow accordingly, and the time derivative of the total 
energy is given by 



d£_ 
'dt 



sp 



(9Tj_ 2 m ^ d q\\ q± d /I m q± 



(73) 



which is the direct correspondence to the hjdg/dt) form in Eq. fl5^ for the delta-f model. 



15| . but now it is a result, not a construc- 



This is the same energy equation as given in Ref. 
tion. 

Following these we can find the moments over h and g defined in Eqs. f l22|) by inserting 
Eq. for 6f and then straightforward evaluation of the velocity space integrals using the 
definitions in Eqs. dnSED]) for the potentials. The moment list for h in entirety is 



nu\\ 



dWzh 



nT\\ = / dW {mz^ -T)h 



nT± + neVtc = J dW (wB -T)h 

gil = I dW {mz^ - 2,T)- h 
J 2 

q^ = j dW {wB-T)zh 



(74) 



similar to Eqs. (143 p as only n and Tj_ are affected. The moment list for g in entirety is 



n = I dW g 

nu\\ + n — Ac = idWzg 
mc J 

nf\\ = j dW {mz^ -T)g 



riT^ = J dW [wB -T)g 

q\\= j dW{mz^-^T)'-g 

q± + riT—XG = f dW {wB - T)z g 
mc J 



(75) 



where due to the extra factor of z it is u\\ and q± which are affected. The moments over g are 
the quantities appearing under d/dt in the gyrofiuid equations, while those over h appear 
under the derivatives in the linear terms. 
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The salient result of this section is the way the interaction Lagrangian and the Hermitian 
property fix the gyroaveraging operations for the potentials once they are decided for the 
moment variables, and vice versa. This relationship is what guarantees an energy conserving 
model using this particular procedure. And given these constraints, the free energy theorem 
follows naturally. 

VI. THE GYROFLUID MOMENT EQUATIONS 

With the above in place, the straightforward derivation of the gyrofluid moment equations 
is essentially determined. The form of 6f in terms of the moment variables, the definitions 
of the latter (hence which moments to take), and the closure rules are already defined. 
We simply take the moment list defined in Eqs. ( H3|) and apply each one separately to 
Eq. (125]) . With one exception, all the steps follow directly. The exception is the parallel 
magnetic nonlinearities, in which both field variables appear, each with a factor of Jo, and 
the velocity space integral must be decided. But as we will see, the form this must have is 
already determined by the requirement of energy conservation. 

A. Toroidal curvature, magnetic divergence effects 

The delta-f form of the gyrokinetic equation splits the curvature and grad-B drifts from 
the rest, as these are purely linear terms. These two drifts are combined, such that in each 
case we take a combined moment over mz'^ + wB. The quantity operated upon is h, not 6f, 
so that the field potential moments also appear. Only is involved, and Jq appears only 
with (f). Moreover, only derivatives over x and y are involved, while B depends only on s. 
These terms are therefore found by simple evaluation of the moments over mz'^ + wB times 
h. These terms conserve energy separately from those arising from other brackets, as in the 
delta-f Vlasov equation (in Eq. |28] h combines with itself under /C as a pure divergence). 
The 4th moment (w'^B^) over Jo0 in Eq. (!62|) is determined by this requirement, as noted 
there. After evaluation of the moments, the curvature terms appear as 




(77) 



(76) 
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^(-^]=--- + ^]c(\+4^] (80) 
dt\e riTj \2 riT J ^ ^ 

(9/1 mq±\ niT fl^ , /qi\ 

The terms with /C operating directly on the moment variables always appear with factors of 
T/e in these scaled units, reflecting the charge separation effect of "diamagnetic" curvature 
terms. Those acting on (p give the ExB compression effects plus FLR corrections. In the 
parallel flux variable equations there is the additional factor of m/e matching the one in the 
inertia terms. These curvature terms form a closed set which conserves energy separately, 



as detailed in Ref. 
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B. Nonlinear ExB advection 

In the nonlinear brackets [(Jo^), {Sf)]xy both and Sf appear under gradients. The 
moment mixing for state variables is given by 



J dW [(Jo0),(5/)] 



xy 



71 



xy 



5G 



T, 



+ 



xy 



xy 



For flux variables it is similar, 

J dWz[{Jo^),{Sf)] 

J dWz{wB-T) [(Jo0),(5/)] 



xy 



xy 



n 



+ 



xy 



xy 



G, q±]^y + ^G, {nTu\\ + 2gj 



xy 



For the pure parallel velocity moments there is no mixing as no factors of wB appear. 



dW {mz^ - T) [( Jo0), (5/)],, = n 0^, T, 



dW-{mz^-T) [(Jo0),(5/)] 



xy 



xy 



xy 



(82) 

(83) 

(84) 
(85) 

(86) 
(87) 



The nonlinearities involving A\\ properly belong to the parallel dynamics. Since the parallel 
gradient ultimately acts on h rather than 5f in the kinetic model, these nonlinearities involve 
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two appearances of the field potentials {A\\ and 0, in the nonlinear parallel electric field). 
In the ExB advection terms these pieces vanish because they all arise from [(Jo0), {Jo4')]xy 
Indeed, if we apply the rules for h to the above six combinations {(pc with n and with T±, 
as in Eg. all the terms quadratic in (p vanish, as they should do. But for [(Jo^y), {Jo4')]xy 
several field terms survive (the lowest order among them being the magnetic flutter effect on 
the electric field in the Ohm's law, in the fluid sense) and the energy theorem's constraints 
are needed to evaluate them. They are left to the discussion on nonlinear parallel dynamics, 
which follows next. 



C. Parallel dynamics and magnetic nonlinear it ies 

In the gyrokinetic nonlinear bracket both (p and Ay appear in the gyrokinetic potential 
(perturbed Hamiltonian) . Since Ay appears with an extra factor of z, we compute these 
terms separately as they involve different moments. The lowest order A\\ terms combine in 
the nonlinear parallel gradient, 

V|| = f (88) 

to which the r2-dependent FLR effects add and mix moments in the same way as with 
(pG and Qg in fhe ExB advection. However, the d/ds terms also involve the dependence 
B = B{s), so it is useful to consider these separately. 



D. Linear parallel dynamics 

The linear terms arise from the [Hq, h]zs bracket in Eq. (!28l) . whose two pieces combine 
parallel streaming and magnetic trapping effects. Here, B does not commute with d/ds so 
there ar,se extra terms proport.oual to BB/Bs, referred to as .nagnetie pu,np,ng terms, the 
vestige of kinetic trapping in the gyrofiuid model There is also the distinction between 
a parallel divergence and a parallel gradient, according to whether 1/B occurs inside or 
outside of d/ds. As with the curvature terms, the moments are over h, not 6f, so the 
same combinations arise as in the curvature terms. With these extra considerations the 
moment calculations are straightforward. Again, these terms conserve energy separately 
and hence can be considered separately (with the magnetic pumping terms independent of 
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the others). After evaluation of the moments, the hnear parallel gradient/divergence and 
magnetic pumping terms appear as 

|5...._.v;| (89) 

|(i^, + -s,,)=..._vs!!*L±I^-M?£±^l^v5iogB (90) 

ot \c e J " ne e " 

dt\e riTj 2 " e ^ ' 

where V|j denotes the linear part of the parallel gradient. In field-aligned Hamada coordi- 
nates it is given by 

d 

SVij = B'— (95) 

with 5* independent of s. The Jacobian does not appear because it is a flux function 
(function of x only). Under strict delta-f ordering both the Jacobian and B'^ are constants. 
In general the Jacobian enters, so that 

where g is the determinant of the metric coefficients (all three components). It is essential 
that ^/gB^ be a flux function, to preserve V • B = 0. 



E. Nonlinear parallel dynamics 

The lowest-order nonlinear terms follow by consistency from the linear ones, as both pieces 
of the nonlinear parallel derivative Vy = — [AcJ^y act together. However, treating 
temperature dynamics with FLR consistency adds FLR nonlinearities to these, involving 
brackets with xg- These involve moment mixing in the same way as for 0g and fla in 
the ExB advection. The moment integrals encountered are basically the same, because the 
integrals over z and w separate, and the parallel dynamics is merely one order higher by the 
factor of z in the moment hierarchy. 
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Terms arising from [( JqAh), {Sf)]xy involving Aq simply follow from the linear ones arising 
from V|j(5/). Those involving xg raise the moment level by one in w in the same way as 
in the ExB advection terms. The only subtlety is the one involving moments over the 
field nonlinearity bracket [Jo^y, ^o^lxj/ which is unique to the nonlinear parallel dynamics 
(the corresponding terms in ExB advection vanish trivially with [Jo0, JQ4>]xy)- These terms 
appear in the flux variable equations, 

d /I m \ r z 

_ {.Aa + -U,)=... + JdWnr. — -IJ,A„ MU (97) 

d fl mq^\ f , wB-T F^^ z ^ ^ 

since in the state variable equations they vanish due to the odd symmetry of / dW z. The 
difficulty is that the velocity space integration is no longer associated with only one of the 
quantities appearing under spatial derivatives in the bracket. The solution is to apply the 
moments to the state variable terms in the {6f) representation in Eq. (jH]) first. These give 
d /I m \ 1 1 ~ 

dt\enTj 2eBo^ ^' "^"^^ ^ ' 

d fl m q±\ 1|. ~, 1 

50 neho 

Now we apply the combination rules for h given in Eq. fl56|) . so that these become 

d / \ TTL ^ \ 1 ^ 

— -Ag + —u\\ =■■■ + — ^[Ag, {p\\ + ne(f)G)]xy 
at \c e J neBo 

+ ^[XG,{f± + enG)U (103) 
e-Do 

dt\e nTj 2eBo^ ^' "^"^^ ^ ' 

d fl m q±\ 1 r . ^ M 

1 _ 2 ~ 

+ — ^[Xg, {p\\ + ne(j)G)]xy + —^[xg, {T± + eQG)]xy (105) 
netjQ e±)o 

Hence we have determined that the rules for moments over two operations by Jq on field 

variables are 

^ J c?WF^[(JoA||), (Jo0)] = [(EiAii), (ri0)] + [(r^Aii), (120)] (106) 
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^ / rfW^F^^[(JoA||), (Jo0)] = [(Fi + T,)A\\,T,<P)] + [r^Aii, (Fi + r^)^] (107) 

It is simple to show that these vanish if Ay — )■ 0, as they should do. 

In the state variable equations the A\\,(f) combinations do not appear, and the nonlinear 
parallel divergences of flux variables are left as 

d fi 1 1 

It = ■ ■ ■ + STk'-""-'^'- ^ 

Here, there is crosstalk between n or Ty and gj_ due to the fact that the corresponding 
wB — T moment does not vanish when a factor of Jo is present. Again, this set of terms 
conserves energy separately from all the others. By contrast to ExB advection, however, the 
conservation is between state and flux variable sets, not just within each of those sets. When 
the flux variables are kept as time dependent variables, the level of moment mixing within 
state and flux variable sets remains consistent (hence velocity with density and heat flux 
with temperature, with one heat flux for each of perpendicular and parallel temperatures). 
With these terms determined, the dissipation free part of the equations, the part which 
conserves free energy exactly, is now closed. Under the strict delta-f ordering, with the 
geometry dependent upon s and with nonlinear derivatives only in {x,y}, the constraint 
V ■ B = is maintained in the nonlinear terms as well. 



VII. DISSIPATIVE EFFECTS 



The addition of dissipation is essentially done by hand, as in most of the moment ap- 
proaches (cf. Ref. l|), as it is not part of the original Lagrangian/Hamiltonian formulation. 
Even at the gyrokinetic level, if collisional effects are very weak then either the computation 
must resolve very thin striations in phase space (cf. Refs. 38N41|). or it must cut them 
off via a hyper diffusion in velocity space which accounts for what is essentially the same 
;hing as a "high Reynolds number" situation viz-a-viz Fokker-Planck collisional diffusion 
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43j. Here we note that the basic mathematical properties of the one-dimensional (ID) 



Vlasov-Maxwell system of the earlier references are the same as the shear- Alfven subset 
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of the electromagnetic gyrokinetic turbulence discussed in the more recent ones. The gen- 
eral gyrokinetic coUisional process is being developed on more firm mathematical grounds 
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45| , but it is mostly simplified models that are currently in use j46[ , when collisions are 
used at all. 

Dissipation effects due to both Landau damping and collisions are added to the gyrofluid 
equations at more or less the same level of sophistication. While it is possible to mimic the 
inear response of the plasma dispersion function in the linear limit in homogeneous geometry 

m 

or for simplified instabilities in toroidal geometry a general fitting approach was 
shown to fail even for linear instabilities close to and away from threshold, by the necessity 
to use a different fitting matrix for each case [l]. Hence, Landau damping was inserted 
into the heat flux equations as a direct-damping model in Ref. ^5)], such that the Alfven 
damping response is adequately captured without impacting any of the conservative transfer 
processes (i.e., by changing their coefficients, as was done in Ref. [2|). 

Collisions are another matter, as there is a well formulated fluid limit to which the 
equations should relax if the collisional frequency becomes large compared to any parallel 
transit or nonlinear advection frequencies — the Braginskii equations jlo|. Typical fluid 
models in the drift freque ncy reg ime keep the parallel dissipation effects (resistivity, thermal 



forces, parallel viscosity) 47H49[. Gyrofluid equations, functioning in the same regime, keep 



these same processes as wellj4, |l5|. The only complication is to keep them consistent 
with the anisotropic temperature model, and to formulate the thermal forces such that they 
disappear naturally when the collision frequency drops to zero j47]. At the drift kinetic 
or gyrokinetic level, the dissipation is covered by the 1-D part of the equations describing 
parallel dynamics, and a simple procedure has been shown for the drift kinetic equation for 



electrons in Ref. 50| which we can use here. This is to set up a Chapman-Enskog expansion 
around a simple Lorentz collision operator, obtain directly all the necessary terms, and then 
change the coefficients such that the result agrees with the Braginskii model. Here, the only 
extension of that was to start from a bi-Maxwellian with given perpendicular and parallel 
temperatures 

= n{2nTii/my^/^{2TTT^/m)'^ exp{-wB /T^_ - mzV2T||) (111) 
and then to solve the correction equation 

27 



for f^^\ where coordinates {v, are given by mv'^ = 2wB + mz'^ and ( = z/v, and ul is the 



Lorentz colhsion parameter. Following Ref. 50|], at lowest order there are no flows (in the 



local rest frame of the ion fluid), so the consistency conditions by which d/dt is eliminated 
are trivial (small parallel force imbalance, no divergences). Also, a flnite A\\ is neglected by 
assuming the resistivity is sufficiently large, and toroidal drifts are neglected by assuming 
R 3> L±, where L± is the profile scale length. Having solved for f^^\ the heat fiux variables 
gilll and q± and, for electrons, the parallel current Jy = —neeUe\\ are evaluated directly. This 
yields collisional formulae for the fiux variables 

ir^^^J\\ = UeCEw + V||J9e|| + aeU^ \\Te\\ (113) 

geiiii + 1.28(ge|||| - 1.5ge±||) + |«e — ^|| = -|fi:e^neV||Te|| (114) 

2 T 2 T 

qe±\\ - 1.28(ge|||| - l-5ge±||) + -ae — J\\ = --Ke ^^eVllTe^ (115) 

5 e 5 niei'e 

(116) 

as well as anisotropy dissipation corrections to the temperature equations 

„^|i,....ii.(r,-Tj (118) 

Here, the numerical coefficients resulting from the Lorentz model which are kept are written 
explicitly, while the coefficients to be substituted with their Braginskii values are written as rj, 
tte, i^e, and TTg, for resistivity, thermoelectric coupling, thermal conduction, and (eventually) 
parallel viscosity, respectively. As the subscripts indicate, the calculation itself is done for 
electrons. For ions, the thermoelectric coupling is set to zero, and then Ki and tTj are given 
their Braginskii values (for several ion species if desired). 

To treat the transcollisional situation the thermoelectric coupling is reformulated (ob- 
viously, it should vanish as Ue — )■ 0). The thermal force itself, ae'n,eV||Te||, is substituted 
using the heat fiux formulae (in this, the anisotropy is assumed to be small for large z/g), as 



detailed in Ref. 
fiux variables, 



47] , so that the electron force imbalance is entirely due to dissipation of the 



/I Ml, \ 



J|| tte / gellll + qe±\\ , -^|| 
V— + —\ 7f^ + "e — 



(119) 
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having The heat flux formulae are treated in a similar fashion, 



PeVllTell = -meUe 



5/2 



2 



5/2 



geiiii + 1.28(ge|||| - l-5qe 
e±\\ - 1.28(ge|||| - l-5ge_ 



N 3 Te ^ 

+ -ae — J\\ 
5 e 

5 e . 



(120) 
(121) 



with gradients on the left and dissipative damping terms on the right. The right hand 
sides of Eqs. f lll7H121"l) represent a set of terms which are added to the right hand sides 
to the corresponding gyrofluid moment equations. The resistivity combination in Eq. flll9p 
subtracts from the right hand side of dAc / dt. The temperature terms in Eqs. flll7|118p add 
to the temperature equations as shown. The heat flux terms in Eqs. fll20fl21l) add to the 
equations for ge|||| and gexy, respectively. The ions are done the same way as the electrons 
except for the resistivity terms (which are the same for all species) , and with the coefficients 
Ki and VTj the appropriate ones for each ion species, with a set to zero for the ions. 

Parallel viscosity does not explicitly appear in the equation for the parallel velocity mo- 
ment but instead results from collisional dissipation of the difference between perpendicular 
and parallel temperatures. This difference is of course one and the same with the parallel- 
parallel component of the viscous tensor 10|. Either one keeps a viscosity term in the 
parallel momentum equation or one keeps track of thermal anisotropy with the dependent 
variables, but not both. There is some discussion of this in Refs. [l, 2|, and it is the ulti- 
mate reason that no Landau closure dissipative term should appear in the parallel velocity 
equation itself jsl. Such a two moment dissipative closure has been given [55, but especially 
for electrons it has the undesirable property of mimicking a resistivity which is a factor of 
Ve/vcQR too large. Hence the Landau closure itself appears only in the heat flux moment 



equations, following Refs. flS^ . and the parallel viscosity is naturally given by the col- 
isional dissipation of thermal anisotropy if the collision frequency is sufficiently dominant 
51 1- The Chapman-Enskog procedure gives the coefficient z/^/S for species as in Ref. jl|. 



The coefficient tt^ then gives the correct viscosity coefficient for species z according to the 
collisional fluid derivation jlo| . 

The correspondence of these equations to low frequency fluid drift equations including 
nonlinear polarisation and collisional dissipation processes in the appropriate ( "Braginskii" ) 
regime has been shown elsewhere and a summary is given in Sec. IVIIIDl below. 
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VIII. RESULTING GYROFLUID EQUATIONS 



The results of the above derivation are essentially the same as those of Ref . 15| . The two 
new elements are the derivation path and the inclusion of direct finite gyroradius effects in 
the part of the nonlinear dynamics representing the fluctuations in the magnetic field. The 
derivation is now firmly grounded within the underlying gyrokinetic theory to the extent 
that the representation of the gyrokinetic distribution function in terms of the gyrofluid 
moment variables is explicit and the energy theorem not only remains intact but has itself 
been used to determine the closure rules. The part of the gyrofluid moment equations which 
arises directly from the gyrokinetic model and is exactly conservative is given by equations 
for the gyrocenter density, 

"^""^ + [fie, T.J + fiV|3 - J = IC ( r.M±^ + <t^G + ^] (122) 



dt ' ' " B -'^L^-'^^-J V 2 2 

parallel velocity, 

/3e-^ + + A^^^G, (lzl\ + V||(0G + rzVz\\) - /3e[XG, (^G + T.T;j_l)] 

+ iS^G + r,r,x - r,r,||)V|| logi? = r,/x,/C (^l^ll±^M±3ll_"^ (123) 



parallel temperature. 



1 ^^-11 + 5V||^^^4^ - + (w.ii + g.±)V|| log 5 



1 dt " B 



2 +T ^'''^ 



perpendicular temperature, 
dT 



dt 



-- + [no, {n, + 2T,x)] 
+ SV||^ - PelxG, {uz\i + 2g,x)] - (m.ii + g.±)V|| log 5 



parallel/parallel heat flux. 



dt 

and perp/parallel heat flux 



f^z^ + ^r.ViiT.ii = r./i./C (^!M±^ ) (126) 



o dxG , dqz^ , ro ^ i o m 
Pe-^ + 1^^—^ + I^zPg, {uz\\ + 2g^±)J 
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+ {Qg + rzT^± - r^T^ II )V II log 5 = r^/x^/C ^^^-^ 



The lowest order nonlinear advective and parallel derivatives are given by 

d d 



V|, =b- V-/3e[A 



G, 



(127) 



(128) 



with b the unit vector of the unperturbed magnetic field. 
The nonlinear brackets are given by 



[/,^7] = V/-Fo-V^ 



5iven V ■ Fq = and Fq : V V = 



(129) 



where Fq is a divergence free, antisymmetric tensor as specified, in the particular normali- 
sation being used. Similarly for the curvature operator. 



/c(/) = rvj 



given Vi/C' = 



(130) 



A typical case is to leave the drift scale ratio Ps/Lj_ out of the normalisation and put it into 
Fq and /C. Conventional normalisation is to fold it into the normalisation, keep to strict 
fiuxtube ordering, so that [/, g] = f,xg,y — f,y9,x in the conventional linearised gyro-Bohm 
version as in Refs. . Both versions are covered in Ref. \\1 
The FLR reduced potentials are given by 



5G 



Tic 



;i31) 



= XG = Ta^ii (132) 

in terms of the field potentials and Ah. The associated field potential equations are given 



by 



Tin, + T2T,± + 



Tn-l 











(133) 
(134) 



for polarisation and 

vl^ii + [riu,|| + r2g, 

z 

for induction. In a computational model the induction equation is actually solved using the 
combinations under the d/dt in Eqs. fll23|127p 



^ (r? + r^)] - vi) All 

= ^ — Ti (/^e^G + ^^zUz\\) + if^eXG + fJ'zQz 
z l^z 



(135) 
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as these quantities on the right hand side are what are actually advanced. 
The constant parameters 



(136) 



ZTe Zm B 

given the background charge density, temperature/charge, and mass/charge ratios, nor- 
malised to electron and deuterium values. The pressures are linearised, so that 



Pz 



Pz 



(137) 



under gradient operators. The profile gradients for all the state variables are included as 
part of the dependent variables; these may be split in a traditional manner without loss of 
generality. 

To these equations are added the coUisional dissipation model, 

dAc du,. 



2 dt 



+ ^J^z 



iz\\ 



dt 



37r. 



iTz\\-T, 



V-Al + — + qe± + CteJl 



dt 

dqz\\ 
" dt 

dqzi 
" dt 



+ -^(T,||-T,^) 
37r-. 



- I^z^z 



5/2 
Kg 
5/2 

Kg 



Qzw + 1.28(g2|| - 1.5g^^) + -a^Jy 

2 

qz± - l-28(gz|| - 1.5q^±) + -otzJ\\ 

5 



with parallel current 



and with numerical coefficients 



^11 = E 



(138) 
(139) 
(140) 
(141) 
(142) 

(143) 



sp 



= 0.71 Ke = 3.2 TTe = 0.73 T] = 0.51 



for electrons and 



ai = 0. 



Ki = 3.9 



vr,- = 0.96 



(144) 



(145) 



for singly charged ions. Note the appearance of r]\\ in the li^y -equation for all species, since 
the resistivity essentially adds to dAc / dt, and that the thermoelectric coupling between Jy 
and the heat fluxes affects only the electrons. For other charge states the corresponding 
coefficients for ions may be found in Ref. |lo| . 
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Finally, the Landau damping model is added 

dqz\\ 



/i2 



dqz 



dt 



■ - liz^^z (l - 0.125g2i?2vj) 

■ - [iz\f^z (l - 0.125g2i?2vj) 



(146) 
(147) 



which gives a finite-difference compatible version of the original by Hammett and Perkins 
and their successors |l, 2|, as explained in Ref. jist . This completes the description of the 
six- moment gyrofiuid model ("GEM" from Refs^jJ, [l4, [xsj), now extended to incorporate 
finite gyroradius effects in the nonlinear magnetic fluctuation dynamics. 



A. Gyrofiuid equations for conisionless reconnection 

Simplified two dimensional models are often used in studies of collisionless reconnection 
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-|55|. In terms of the geometry only the dynamics perpendicular to a prescribed 
guide field is retained; in the language of this work this means only the nonlinear brackets 
are kept, with the linear parallel derivative incorporated into the magnetic nonlinearities, 
assuming a Cartesian coordinate system which is not aligned to the component of the mag- 
netic field described by the shear. The standard "four field model" is an isothermal two 
fluid version of the equations keeping parallel velocities and densities. In fluid language it 
is 



30| 



drie 

"dt 

du\\ 



V||J|| (148) 
V||(J||-«ll) (149) 



dt dt 
where the vorticity and current are given by 



dAw dJw 

/^e^+/^e-3r = V||(ne-0) (151) 



= V10 Jii = -VIA,, (152) 



and the nonlinear derivatives are given by 

d d , 

dt = m^^^^^ 



] (153) 
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This is obviously related to the above gyrofluid model, with cold ions and with isother- 
mal electrons, and no FLR effects. The simplified "two-field model" is then prescribed by 
neglecting u\\ and setting ng = Q. The four field model is equivalent to the two-moment 



simplification of the gyrofluid model, with Jy 

drii 
It 

^ dA\\ duii 



/3e 



dt 



dt 

(if II 

'~dt 



u\\ — f II and = rig — rii, so that 



-V||f|| 

-V||M|| 

-V||0 
-Viiid - rie. 



(154) 

(155) 

(156) 

(157) 
(158) 



up to some /ie/ /ij corrections. Putting the FLR effects back in (with constant background 
temperatures) we then have the two moment equations for each species, 

dn. 



dt 

dAc , du,.\\ 



dt 



dt 



-VyU^II 

-V||(0G + r^n^ 



(159) 

(160) 
(161) 



with nonlinear derivatives, 



d _ d 
di~ di 



+ 



V|| = -f3e[AG, 



and FLR potentials. 



and with polarisation, 



Tie 



A 



sp 



G 



Fn- 1 



Ti^ll 



and induction (cf. Eqs. I134fl35p . 



e p2 
1 



- Vl Uii = E — Ti (PeAa + fi.u, 

J z l^z 



(162) 



(163) 



(164) 



(165) 



Hence, we find that the four field model has this obvious generalised FLR version, through 
the two gyrofluid moments, and the two-field model is basically the same thing with the 
ion gyrofluid moment variables neglected on the basis that self consistent gradient drive is 
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absent and Pe/ ^J^i ^ 1- The foregoing correspondence was shown in Refs. ll5|, in the 
context of tokamak microturbulence. 

The obvious next step for reconnection with finite electron gyroradius modifications to 
the nonhnear magnetic field dynamics is to restore all six gyrofluid moment variables, 

' + [0G, n,] + [fie, T,^] = l3e[AG, + /3e[xG, Qza] (166) 



dt 



OAg duz\\ r . 1 !■„ 



l3e[AG, (0G + r.p,|i)] + /3e[XG, {^G + T.T.x)] (167) 



l9T,ii 1 



2 dt 2 

dt 



+ o[0G,T,|,] = /3e[AG, (W.ll + g.ll)] + /3e[XG,g.J (168) 



+ [0G, r,^] + [fie, (n. + 2T,^)] = (3,[Ag, + /3e[XG, (m.|| + 2g.^)] (169) 

/^^^^ + = ^1"J^G,^^||] (170) 

= r,/3e[AG, T,^] + /3e[XG, (0G + r.p,||)] + /3e[XG, 2(1)g + TzTz^)] (171) 

with FLR potentials and polarisation and induction equations as in the full six-moment 
gyrofluid model (Eqs. I131H134I) above. 



B. Constraints on magnetic geometry 

At this level the derivation of a set of equations actually used in a computation depends on 
assumptions made about the geometry. For 2D reconnection studies in a model containing a 
guide field the typical case is that a background homogeneous magnetic field Bq is prescribed, 
the coordinates of the computational domain describe the plane perpendicular to Bq, and the 
components of the magnetic field in this plane, labelled xy, are determined by the evolving 
potential A\\. For 3D magnetised plasma turbulence studies the coordinates are usually 
aligned to a background magnetic field with both curvature and shear. The derivation of 
gyrofluid equations does not depend on the use of either, or other, sort of model. However, 
the most important point is this: energetic consistency must be maintained. In this context 
the only requirement is that the tensor Fq defining the nonlinear brackets is antisymmetric 
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and divergence free, and the vector /C* defining the curvature operator is divergence free. 
Although the equations were derived using strict fiuxtube ordering, that was made necessary 
by the dependence of F^'^ upon B and therefore its spatial dependence. That was simply an 
energetic consistency constraint. Once the equations are derived, one need only maintain 
that same level of energetic consistency, and this is done by retaining the required properties 



for Fo and /C. This makes a globa. 



geometry is given in Refs. 
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geometry model possible. Local and global field aligned 



561 ] and the necessary constraints of global consistency of 



luxtubes and the transforms necessary to obtain global mode structure are given in Refs. 



36l | and |37|, respectively. An example of a global tokamak geometry model which relaxes 



the ordering on the coordinate derivatives, while keeping it on the equations themselves, is 
to recast the curvature operator in terms of a bracket. 



/C(/) = 2[logi?,/] 



(172) 



with log R a scalar function of the coordinates, and the brackets given as in Eq. 11291 with 
log-R and the components of Fq, along with the metric and magnetic field component i?^, 
describing the geometry. A version of this is given in Ref. 37| and is suitable for global 
computation of the tokamak core, or the edge region within one scale-length variation of the 
background parameters dependent upon densities and temperatures. The ordering behind 
the equations essentially limits validity to a domain comprising one such set of scale lengths. 
Otherwise, the parameter set ceases to be representative (cf. the discussion in Refs. fl^ . 



57(1, or 58|). These steps are what extend the GEM model of Refs 



model described in Refs. 



51 



59] 



.0 



14 



15| to the GEMR 



C. Global computation and the issue of stratification 



Global computations with this delta-f gyrofluid model and its predecessor in Ref. [1 
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m and 



are not only possible but have been underway for some time, both for edge 
for core [59| cases. With the background current profile added to the electron parallel 
velocity dependent variable Ue\\ — ?■ "UeH — Jo/i'^e^), where Jq = Jo{x), with UeC a normalising 
parameter, the model is being used for studies of the self consistent interaction between 
dominantly ion temperature gradient driven electromagnetic core turbulence with magnetic 



islands 



6l| . One often thinks of fiuxtube computations as being defined by periodic radial 
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boundary conditions. In fact they are defined by taking the radially local approximation 
on the flux surface geometry, so that all metric quantities depend on s only except for the 
appearance of shear, and then assuming that d/dx, d/dy ^ d/ds in the derivatives, (cf. Sec. 



IB). 



14 



'luxtu 
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DC computations without periodic boundary conditions are common {e.g., Refs. 



58 



62 



64|). Although the delta-f gyrofiuid equations are initially derived using 



fiuxtube approximations on the gyrokinetic model , in the resulting gyrofiuid model global 
geometry may be restored as long as the model drift tensor used to define the brackets is 
divergence free (cf. the previous subsection; this step cannot be taken directly on the delta-f 
gyrokinetic model, due to the dependence of F^^ upon B and the need to commute F*^ past 
derivatives to conserve free energy). 

The other issue faced by truly global computations is stratification. Although delta-f 
models can use global geometry and do global simulations, what they cannot represent is 
a change in physical parameters across the domain. Again, this is for energy conservation 
reasons: the physical parameters in the delta-f equations (including those in F^^) must 
also commute with the nonlinear bracket operations, including the /C terms. A stratified 
nonlinear model will require that all locally varying parameters {i.e., temperature, beta) be 
determined (three-dimensionally) by the dependent variables. Then, the conserved quantity 
is no longer free energy but the Noether energy, as explained in Sec. Ill CI For a gyrofiuid 
model to be able to represent stratification, analogous to nonlocal fluid Braginskii equations 
of the sort given in Ref. 65|, it must either be derived from a total-f gyrokinetic model such as 



given in Sec. [TTl or the simplified version without FLR effects usee 



;o treat equilibrium flows 



19| . This is still ongoing 



variables and the simultaneous 
18, 13, which is the closest at 



in Ref. [66|, or directly from a model Lagrangian as in Refs. 18 
work as the need to have parallel heat fluxes as dynamica 
incorporation of FLR effects requires extension of Refs. 
present towards a complete total-f gyrofiuid model. 

With such a stratified model one can anticipate to capture the spatial transition between 
regions of different physical character and perhaps also temporal transitions. Of course, to 
treat such phenomena one requires a solid foundation within energetic consistency. 
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D. Relation to dissipative fluid drift equations and Reduced MHD 



A detailed analysis of these dissipation terms under the Braginskii limits (small gyro- 
radius, strong collisionality, and an implicit assumption that specific heat fluxes are small 
compared to fluid velocities) showed that they recover the forms used in reduced (low fre- 
quency, neglecting compressional Alfven dynamics and expressing perpendicular flows as 
drifts rather than directly as dependent variables) Braginskii equations. Polarisation heat 
flux effects due to finite perpendicular inertia were also recovered. By extension, the re- 
duced magnetohydro dynamic (MHD) equations are already well known to be a subset of 
the reduced Braginskii equations. These results were obtained in a previous work [67|, and 
are summarised here. 

The Braginskii regime is defined by the assumptions used by Braginskii to obtain his col- 
lisional fluid equations 10| , whose low frequency regime still represents the most commonly 
used model for nonlinear edge turbulence computations. The assumptions are: (1) domi- 
nance by collisions, such that both ^ d/dt and z/^ ^ V^^Vy, where z/^ is the like particle 
collision frequency and Vz is the thermal velocity, (2) long wavelength, such that p^V^ ^ 1, 
where is the thermal gyroradius, and, much less commonly understood, (3) small specific 
heat flux, such that ^ Pz^z-, implicitly assumed in the use of a drifting Maxwellian to 
lowest order and obtaining q only through first order corrections. In these inequalities the 
subscript z denotes the species, i.e., the inequalities should hold for all species. The most 
important thing to know about tokamak edge turbulence (in particular) is that all three of 
these assumptions are violated, and the problems are especially severe in the ions, not the 



electrons which received most of the early attention [50| . Specifically, z/j is one to two orders 
of magnitude slower than the turbulence, for the longer and shorter wavelength component, 
respectively, and the turbulence vorticity spectrum always extends past A;_lPs = 1 for drift 
wave turbulence (in a plasma with Tj ~ T^. this means a component with k^pi ~ 1 must 
be faced), and especially for temperature gradient driven turbulence, not only are specific 
parallel heat fluxes larger than velocities in fluctuations, but the diamagnetic compressional 
effects in tokamak geometry scale with the gradients and therefore the heat flux effects are 
stronger in both the diamagnetic and polarisation effects in the ions. The latter indicates 
extending the fluid models to treat heat fluxes as well as velocities in the stress tensor effects 
that represent nonlinear polarisation, as introduced and discussed in Ref. 68|. 
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The gyrofluid equations represent all these extension effects as well as the basic Bragin- 
skii ones automatically. Brieffy, nonlinear finite gyroradius effects in the gyrofluid model 
reduce term by term to the nonlinear polarisation velocity and heat flux effects in an ex- 
pansion keeping lowest order fc^p^ effects in the ExB advection and in the equivalence of 
representation implied by the polarisation equation. For an isothermal model, one can start 
with the definition of 0g, solve the polarisation equation for rii in terms of ne, expand these 
expressions one order in k^pl and insert them into the nonlinear advection derivative for rij, 
and do some operations involving converting V^^ inside brackets into general divergences. 
The result is the fluid density equation with the nonlinear polarisation terms. With the 
temperatures the procedure is more complicated and involves manipulation of the VLq terms 
and derivation of the space/gyrofluid temperature representation from the moment hierar- 
chy. This procedure recovers the temperature gradient parts of the polarisation velocity and 
also the polarisation heat flux terms. 

The other correspondences involve collisional dissipation formulae and are simpler. The 
parallel viscosity is proportional to the temperature anisotropy 5Tz = Tz\\ — Tz±. The two 
temperature equations are used to form a time dependent equation for 6Tz. The velocity 
divergences (Vy and /C terms) and collisional dissipation are put on the right side, and the 
time derivatives, nonlinearities, and heat flux divergences are put on the left side. The 
Braginskii regime consists of assuming the left side to be small, recovering the collisional 
formula for parallel viscosity (since we set the vr^ coefficients accordingly). As noted, this is 
severely violated for ions in edge turbulence. This result was also obtained in Ref. 51 1. 

The parallel heat flux in a collisional model is the sum of the perp-parallel and parallel- 
parallel components. The equation for these are added, the temperature gradients and 
collisional damping terms are put on the right side and all other terms are put on the left 
side. The Braginskii regime consists of assuming the left side to be small, recovering the 
collisional formula for parallel heat flux (for this we set the Kz coefficients accordingly). This 
is also violated for ions, with the ExB advection nonlinearities much larger than any of the 
dissipation terms for typical parameters. This result was also obtained in Ref. 51|. 

Once the Braginskii regime is (formally) recovered, the next steps to reduced MHD are 
well known (this is a matter of the difference between two-fluid models, e.g., Wakatani- 



Hasegawa 



23 



69|, and the one- fluid versions of reduced MHD 70|, LDJ). At the isothermal 



level, leading to what are called four-field models, this was done in detail with both models 



39 



and results for edge turbulence in Ref. 



72| . The basic assumptions are that e(/)/Te is larger 



than any relative fluctuation variable in densities or temperatures, and that compressibility 
effects in Jy are neglected (energetically, the second assumption follows from the first). This 
can be called MHD ordering (pressure gradients are neglected in favour of electric field 
components). 

At the philosophical level, the gyrofiuid model treats polarisation densities instead of 
velocities, but the simplest road back to reduced MHD is to apply d/dt to the polarisation 
equation, apply MHD ordering, and follow all the consequences. Similarly, the Ohm's law 
amounts to neglecting electron inertia effects and applying MHD ordering to the electron 
parallel velocity equation. Reduced MHD without polarisation effects is particularly easy 
to understand (neglect all FLR effects, subtract the continuity equations, and replace the 
gyrocenter charge density with the ExB vorticity using the polarisation effects). A more 
introductory version of this correspondence between models was given in Ref. jl4]. Basically, 
the gyrofiuid model solves the same problems as the reduced Braginskii fluid model, without 
the more damaging assumptions of the latter. 



IX. SUMMARY AND COMMENT 



The gyrofiuid model corresponding to what is known as the "delta-f" gyrokinetic theory 
now has a derivation path which starts from first principles. It depends on the ordering 
used to obtain the delta-f forms of both the gyrokinetic equation and the associated field 
potential equations. It further depends on a list of moment variables, the description of the 
underlying distribution function in terms of these variables, and essentially one assumption 
involving the moment closure of finite gyroradius (FLR) corrections. The energy theorem 
descends from both the total- and delta-f gyrokinetic versions and is used to fix the rest 
of the undetermined quantities in the gyrofiuid model. The resulting gyrofiuid model is 
now consistent at the level of the best delta-f gyrokinetic models, especially in terms of an 
energy theorem which is properly conservative in all the reactive effects (compressibility, 
coupling to flows and MHD, etc.), and also in terms of the responses of heat fluxes to 
temperature gradients being carried at the same level of sophistication as the responses 
of flows and currents to density gradients and electric fields. For dynamics of magnetised 
plasmas driven principally by temperature gradients this is the minimal fluid model. For 
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dynamics more generically pressure or current driven the isothermal version 1^, easily and 
consistently obtained by setting the temperature and heat flux moment variables, the second 
FLR operation and the thermoelectric collisional effect (through ag) to zero, becomes 
the minimal model. Hence, the decision of what constitutes the minimal model depends on 
the problem being considered. However, the procedure given herein is useful generally in 
constraining the derivation of any gyrofluid model based upon, as it must be, an underlying 
gyrokinetic model which is energetically consistent. 
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